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ABSTRACT 

We study the implementation of Monte Carlo renormalization group (MCRG) in mo- 
mentum space. This technique is most efficient when used in combination with a Fourier 
accelerated Langevin algorithm. As a benchmark we calculate the critical exponents v 
and r\ in the vicinity of both the gaussian and the Wilson fixed point in A0|. The re- 
sults are very competitive with alternative analytical methods and require a moderate 
computational effort only. 
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1. Introduction 

Montecarlo renormalization group [1] is a powerful technique to determine critical expo- 
nents in numerical simulations. In this paper we implement this technique in momentum 
space (as opposed to real space) in conjunction with a Fourier accelerated Langevin al- 
gorithm, used to generate the set of field configurations upon which the renormalization 
group transformation is performed. Our primary interest is its application to the crumpling 
transition in random surface models [2] but we have found the results for a Z2-symmetric, 
one-component, scalar field theory in 3 dimensions interesting in their own right since most 
MCRG studies have been performed with spin models. These results are reported here. 

Fourier accelerated algorithms were introduced by Batrouni et al. [3]. The method 
uses a stochastic Langevin equation in momentum space, partially eliminating critical slow- 
ing down. Being already in momentum space it is natural to define a Kadanoff transforma- 
tion there, discarding the high frequencies of the system and keeping the low frequencies, 
which are the relevant ones near criticality. The procedure is close in spirit to Wilson's 
original proposal[4] and it apparently avoids some difficulties that arise with some 'un- 
physical' renormalization group transformations in real space. The method we discuss also 
possesses the desirable property of being exact in the free case. And, last but not least, 
the present method delivers very accurate results for the critical exponents v and r\ at the 
Wilson fixed point with little computational effort. Sub leading exponents at the Wilson 
fixed point are of course harder to get and we do not have accurate values for them yet, 
but they do not seem to be hopelessly beyond reach. 

2. Numerical Renormalization Group. 

Just to fix the notation, let us briefly review the method. Our starting point is the bare 
action 



If we perform a renormalization group transformation with a scale t we end up with a 
similar action, but with renormalized couplings and operators. The new couplings will be 
highly nonlinear functions of the old ones 




(2.1) 



a 



(2.2) 



Near the fixed point we can linearize such function 




(2.3) 



The largest eigenvalue Xh of T a/ 3 will give the critical exponent v 

v = hit/ In Xh 
To compute T a p we use the chain rule 



(2.4) 




(2.5) 
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The expectation value of a renormalized operator can be computed either directly with 
the renormalized action or, using the renormalization group transformation to express the 
renormalized operator in terms of the bare fields, with the bare action. Thus 



(O^O b )-(O^)(O b p ) (2.6) 



= (0;0$ - (OWOft (2.7) 



All one needs is a set of statistically independent field configurations to determine these 
correlators and, in turn, T a p. MCRG has been most successful when applied to Ising spins 
in 2 and 3 dimensions [1,5]. See [6] for other applications. 

The Kadanoff transformation is largely arbitrary. We shall consider a transformation 
with t = 2, reducing the volume of our system from N = m x n<i x ... x na points to 
N s = ni/2 x 77-2/2 x ... x rid/2. The renormalization group transformation we will use 
to achieve this volume reduction will be a decimation in momentum space. At each step 
we will simply discard the high p modes. This transformation is very easy to implement 
numerically and it has many advantages, as it will be clear in what follows. 

A renormalization group transformation is not just a Kadanoff transformation. It also 
involves a field rescaling because redefinitions of the fields entering in the action make the 
function (2.2) ambiguous. We must pick an operator (usually the kinetic term) to have a 
definite normalization which will be preserved. Let us work out this rescaling. For a free 
theory the action is 

S=l 5>(s) - 4>(x + A)) 2 + ^ £ <P(xf (2.8) 

We set the lattice spacing equal to 1 and assume a hypercubic lattice in d dimensions. 
Define 

^ = Tr E exp(-i$5)0(p) (2.9) 



N 

p 



then 

S = 



2^Z>("P) (^-2£cos|v + m 2 W) (2.10) 
p \ n=i J 



We now discard the modes with > n/2 and re-express the action in real space via 
(2.9), but on the small lattice of size N s = A^/2 d . The result is [4] 

s = lw 5>(*) - + a)) 2 + ^§ E *(*) 2 ( 2 - n ) 

x,/j, x 

To get the same kinetic term we need to redefine <f>(x) = C4>(x) where Q = l/2 ( - d+2 ^ 2 . 
Then the only other renormalized coupling — the mass — is m 2 . = Am 2 , as expected. 
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In an interacting theory we cannot determine the rescaling ( exactly. We shall write 
C = l/2( d + 2 -^)/ 2 where r\ is the deviation in the scaling of the kinetic term with respect to 
the free case. Provided we are close enough to a fixed point, rj is the anomalous dimension 
of the field (j). 

We shall determine ( (and, consequently, rf) selfconsistently by demanding that ob- 
servables computed with the bare and renormalized actions should coincide. The results 
are quite sensitive to the value of r\. Departures from the correct value will quickly lead to 
unacceptable values for the eigenvalues of the T matrix (very large, complex, ...). 

3. The Fourier Accelerated Langevin Algorithm (FALA). 

One alternative to (partially) beat critical slowing down in Monte Carlo methods is 
to use the Langevin equation to construct a random walk in configuration space 



(j)(x,t n+l ) = (j)(x,t n ) + A(j)(x,t n ) (3.1) 
6S 

54>{y,t n ) 



s s 

Acf>(x) = -At J2 efo V) tt? 7T + y/Airjix, t n ) (3.2) 



where rj(x, t) is a gaussian noise satisfying 

(r 1 (x 1 t)r 1 (y 1 t')) = 28 tit ,e(x 1 y) (3.3) 
In momentum space the Langevin equation reads 

5 S 

<f>{p, t n+1 ) = <f>{p, t n ) - Ate{p)F — + y/Ate(p)ri(p, t n ) (3.4) 

and 

(v(Pi,t)v(P2, t')) = 26 tl t'Sp 1 , -p 2 (3.5) 

Given infinite computer resources e(p) is irrelevant. Any will do. However, to beat critical 
slowing down a physically motivated choice is called for. For the free theory, using e(p) = 
D(p, m) (D is the free lattice propagator) is the optimal choice — all modes are updated 
at the same rate [3]. We expect that in the interacting case a good choice will be e(p) = 
D(p, l/£) , where £ is the correlation length in units of the lattice spacing. 

One point to worry about is that for finite At the equilibrium distribution is not the 
one corresponding to S but rather to 

Se , = S + I E A^H^g^ - + CKMf (3.0) 



Our action S will just be (2.8) with the additional term 

^5>(*) 4 (3-7) 

X 
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The finiteness of the Langevin step induces a systematic error of O (At) in S eq . For instance, 
the bare mass gets renormalized by an amount that is approximately linear in At. We 
would therefore be inclined to use a value for At as small as possible, but this requires 
longer runs to reduce the statistical error. 

The latter can be estimated as ^/2r/Ni, Ni being the total number of Langevin 
iterations, fixed by budgetary constraints, and r being the autocorrelation time, fixed by 
At, but depending on the model dynamics too. For the free theory r ~ 1/At (the actual 
proportionality constant depends on the normalization of e(p)). Away from the critical 
point this relation will hold for the interacting theory too; the autocorrelation time will 
be proportional to (At) -1 and approximately A independent. However, as we approach 
criticality the proportionality constant in the autocorrelation time will grow as £ 2 (£ being 
the correlation length). For local Monte Carlo algorithms z = 2, while FALA are 'optimal' 
having z = for multicomponent scalar theories. In our case, one-component scalar A</> 4 , 
the presence of a discrete move 4> — > —(f) degrades the performance to z ~ 1 — 2 [7]. The 
behaviour of the autocorrelation time is shown in Fig. 1. The difficulty with the discrete 
move will be absent in most applications. 

Adding the two errors in quadrature we get that the optimal Langevin time step is 

— 1/3 

proportional to N i . Trial and error also leads to Langevin steps in rough agreement 
with this scaling. A possibility to decrease the systematic error is to use a second order 
Langevin algorithm effectively removing the O(At) corrections from S eq . The resulting 
code is of course slower, poorer statistics partially making up for any gains in precision. 
Our conclusion after some preliminary runs is that it hardly pays to use a second order 
algorithm. 

4. The Gaussian Fixed Point. 

This corresponds to the point A = 0, m 2 = 0. Since we set A = from the start, the results 
we will present are trivial to obtain analytically, but they are a test that any MCRG method 
should pass. Autocorrelations are reduced completely (i.e. r ~ 1) if At ~ 5 x 10~ 2 (even 
though such big Langevin step induces a renormalization on the fields and the mass itself). 
We have checked that Wick's theorem is satisfied with very good accuracy. 10 bins with 
5 x 10 3 configurations per bin were enough to get accurate results. We have used a lattice 
of volume 32 3 reducing it to 16 3 . The rule-of-thumb for At previously discussed suggests 
At = 7 x 10~ 4 which is the value we have used. The operators we have considered are 
O n (4>) = ^ x 4>{ x ) 2ri - an d their respective dispersions a(O n ). No 'derivative' operators 
were considered at this stage. 

A first application of the algorithm taking m 2 = 10 -3 gives the results shown in 
Table 1. We can extract the critical exponent v = .5027(30) (exact v = .5). The second 
eigenvalue is not very accurate. If we are very near a fixed point the renormalized and 
bare action are the same so averages of observables should coincide, that is (0„) = (OQ. 
In this way, by fitting some observables, we can get the anomalous dimension r\. In our 
case, including an anomalous dimension worsens the agreement, so we take rj = 0. 

If we gradually approach the critical point by letting m 2 — > the observables match 
and the subleading critical exponents agree with their exact values (Tables 3 and 4). Our 
best result for the critical exponent is v = .5004(5). Going to much smaller masses is 
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problematic due to finite size effects. The results are remarkably stable with respect to 
the addition or removal of higher dimensional operators. The lesson we extract from this 
simple example is that it it is easy to get the first eigenvalue right, but to get the subleading 
critical exponents one must be very close to the fixed point. 

The errors quoted in the critical exponents are overwhelmingly dominated by system- 
atic effects due to the truncation of the bare and renormalized actions. Statistical errors 
are negligibly small in all cases presented. Errors due to the fact that S eq ^ S are naively 
expected to cancel out in universal quantities. However, regarding this point one must 
keep in mind that, because a finite Langevin step modifies in practice the action of the 
system, S eq is no longer at the gaussian fixed point. Then matching between bare and 
renormalized expectation values is no longer guaranteed (of course, since At is very small, 
the matching is excellent for m 2 = 10~ 5 , but not perfect). For this reason, including many 
operators in the renormalized action does not automatically improve the results; it may 
actually slightly worsen them because small errors in the highest dimensional operators 
eventually affect the precision of the largest eigenvalues. Finally, errors due to finite size 
effects do not seem to be relevant for these values of m 2 . 



Table f : m 2 = fCT 3 , A = 


Operators 


f eigen 


2 eigen 


f 


3.993 




2 


3.977 


1.688 


3 


3.972 


1.775 


4 


3.964 


1.731 


5 


3.967 


1.731 


6 


3.967 


1.731 


7 


3.964 


1.784 


8 


3.967 


complex 


exact 


4.000 


2.000 



Table 2: m 2 = 10" 3 , A = 


O n /N 


big lat. 


small lat. 


1 


0.259 


0.219 


2 


0.202 


0.144 


3 


0.262 


0.157 


4 


0.475 


0.239 


5 


1.104 


0.465 


6 


3.141 


1.103 


7 


10.59 


3.077 


8 


41.47 


9.800 
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Table 3: m 2 = 10" 5 , A = 


Op. 


1 eig 


2 eig 


3 eig 


4 eig 


1 


3.999 








2 


3.997 


2.003 






3 


3.994 


2.005 


1.003 




4 


3.994 


1.985 


1.003 


0.402 


5 


3.993 


1.999 


0.999 


0.551 


6 


3.994 


1.989 


1.019 


0.553 


7 


3.994 


1.998 


0.964 


0.584 


8 


3.994 


2.002 


0.930 


0.620 


exact 


4.000 


2.000 


1.000 


0.500 



Table 4: m 2 = 10" 5 , A = 


OJN 


big lat. 


small lat. 


1 


0.304 


0.308 


2 


0.275 


0.274 


3 


0.410 


0.400 


4 


0.850 


0.775 


5 


2.242 


1.935 


6 


7.150 


5.800 


7 


26.65 


20.70 


8 


113.8 


88.27 



5. The Wilson Fixed Point 

The critical exponents in the vicinity of this fixed point are not known exactly. There are 
many different approaches giving reasonably precise estimates of the exponents, e expansion 
[8], High T c [9], exact Renormalization group treatments [10] and MCRG analysis [5,11]. 
The evidence accumulated suggest that no other fixed point exists and give an average 
estimate for the exponents 

v = 0.627 - 0.631 77 = 0.03 -0.04 (5.1) 

As we have discussed in section 3 even using a FALA does not completely suppress 
critical slowing down all the way to z = in the particular example we are considering. 
To guarantee sound statistics we accumulated O(10 5 ) configurations in bins of 10 4 , until 
observables reach a plateau. We also measure autocorrelation times. Since the propor- 
tionality constant relating r and At can be large near the fixed point, the Langevin step 
has to be relatively large. The compromise between statistical and systematic errors leads 
to At ~ 10~ 3 . The observables get some small renormalization but universal quantitities 
should be relatively unaffected. 
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We must fine tune the parameters m 2 < and A so that we are near the fixed point. 
For each A > there is a value m\ < corresponding to the critical surface. Applying 
renormalization group transformations there would drive us to the fixed point. Unfortu- 
nately we cannot do that in the computer since near the critical surface the correlation 
length starts to diverge and finite size effects enter into the game. What we did is to 
perform numerical simulations for a fixed value of A and then sweep over m 2 to find the 
critical surface. For this purpose we simulated sistems with 16 3 , 16 2 x 8, 16 x 8 2 and 
8 3 sites. (We were forced to study asymmetric lattices because fast Fourier transforms 
require the lattice sizes to be a power of 2.) In addition to measuring observables such 
as energy density, specific heat, magnetization, etc. we also determined the correlation 
length by fitting the long distance part of our measured Green functions to a propagator 
of the type l/A + £ -2 . Once the critical surface is located, a MCRG is performed at a 
value of m 2 in its vicinity, but where the correlation length is not too large for our lattice. 
MCRG is done on a 32 3 lattice so the renormalized lattice becomes 16 3 . We repeat this 
procedure for several A values until we observe a near perfect match between the bare and 
renormalized actions. This will be the approximate location of the fixed point. Of course 
the true fixed point lives in an infinite-dimensional parameter space and we are simply 
getting as close to it as possible in the (A, m 2 ) plane. Because the agreement between bare 
and renormalized expectation values at our approximate fixed point is good we conclude 
that the particular renormalization group transformation we are using has a fixed point 
not far from the (A, m 2 ) plane. 

To begin with, we consider only local operators. We study values of A ranging from 
0.1 to 20. For A = 0.1 the first eigenvalue of the T matrix is 3.969, almost the one 
corresponding to the gaussian fixed point. For such range of parameters we are still in the 
domain of attraction of the gaussian fixed point. Then we increase A gradually repeating 
the same stratregy. The largest eigenvalue gradually changes. By the time we reach the 
region 5 < A < 8 observables in the big and small lattices match with an accuracy similar 
to the one around the gaussian fixed point. In fact, the best agreement is obtained around 
A = 6.5 which we take as approximate location of the new fixed point. The agreement 
improves if we allow for a small anomalous dimension. Fitting the 5 lowest dimensional 
operators gives 



The agreement between big and small lattices is illustrated in Table 5 for the case A = 6 
and Table 6 shows the results for the critical exponent v. Our best estimate is 



The error corresponds to the error in the eigenvalue. Notice that also in this case the 
agreement between expectation values deteriorates for the highest dimensional operators. 
This is again a reflection that our bare action (which contains only three operators) is not 
quite at the fixed point, just close to it. Both v and rj are in perfect agreement with the 
ones obtained in [5,8-11]. If we keep increasing A we are still in the domain of attraction 
of the wilsonian fixed point, but we are unable to match observables in the big and small 
lattices, even by fitting an anomalous dimension. Fixing r\ = 0.04 we get v = 0.630(5) for 



r] = 0.03 ±0.005 



(5.2) 



v = 0.628(5) 



(5.3) 
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A = 10., still in agreement with the value obtained previously. Table 7 shows the big/small 
lattice comparison for A = 10. with r\ = 0. 



Table 5: m 2 = -3.32, A = 6 


O n /N 


big lat. 


small lat. 


1 


0.272 


0.267 


2 


0.161 


0.160 


3 


0.130 


0.133 


4 


0.130 


0.140 


5 


0.151 


0.173 


6 


0.187 


0.240 


7 


0.260 


0.380 


8 


0.402 


0.650 



Table 6: v for different A 


A 


V 


error 


5.0 


0.623 


0.006 


6.0 


0.628 


0.007 


6.5 


0.628 


0.006 


7.0 


0.628 


0.005 


7.5 


0.629 


0.005 



Table 7: m 2 = -4.94, A = 10. rj = 0. 


O n /N 


big lat. 


small lat. 


1 


0.238 


0.185 


2 


0.120 


0.087 


3 


0.080 


0.061 


4 


0.067 


0.055 


5 


0.063 


0.059 


6 


0.066 


0.073 


7 


0.076 


0.101 


8 


0.093 


0.153 



We also explored the range of A values up to 20., but the matching gets much worse 
and we eventually hit convergence problems of the Langevin algorithm. It is also important 
to stress that no evidence for another fixed point was found, and even for the high values 
of A available the results are in the domain of attraction of the wilsonian fixed point. 

The method delivers v (and rf) with good precision. Unfortunately, the second eigen- 
value at A = 6.5 is 0.3, which corresponds to u = 1.7, far from the accepted result [10], 
but these are subleading effects. One should really get much closer to the fixed point to 
be sensitive to these. 
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It is very important to check the stability of our results under the inclusion of 'deriva- 
tive' operators. Exact renormalization group treatments did not include them in their first 
studies[10a], but some recent ones do indeed include some[10b]. We can consider their 
effect by including operators of the kind 



We carefully studied the inclusion of these operators up to n = 6 in the computation of the 
T matrix in the range 5 < A < 8, but the first eigenvalue did not change at all. However 
the second eigenvalue changed significantly to A = .82 (that is u> = .28), now undershooting 
the generally accepted value which is u> = .5 — .7. Again the situation in the gaussian fixed 
point repeats itself, it is relatively easy to get the first eigenvalue right but the second one 
needs being very close to the corresponding fixed point. 

6. Final Remarks. 

We single out three avenues along which the method could be improved. 

(1) Iteration of the renormalization group: This would require using a bigger lattice. 
On, say, a 64 3 lattice one could perform two successive renormalization group transforma- 
tions. This would surely improve the matching amongst observables further. In principle, 
if our initial bare action is very close to the critical surface, iteration would bring us closer 
to the true fixed point. In practice, however, we cannot be too close to the critical surface 
due to finite size effects and while we certainly expect an improvement on our results after 
iteration this will not be a dramatic one, probably. In fact, we find remarkable that our 
approach yields such accurate results with so little effort. 

(2) Determination of the renormalization group flow trajectories: The method as 
it stands so far is not accurate enough to determine the renormalized couplings. One 
possibility would be to use the method suggested in [12]. Another one would require a 
combination of MCRG techniques with a Ferrenberg-Swendsen[13] algorithm. The idea is 
to use the Ferrenberg-Swendsen method to perform a numerical simulation in the small 
lattice at a value of the parameters where already accurate results are found, storing as 
many operators as possible. Now apply MCRG paying special attention to the observables 
that have been stored. Then we can compute the observables for values of the couplings 
different from the actual one used in the simulation, and choose the couplings that best fit 
the observables. The couplings define the renormalized action. 

(3) Improving the accuracy of subleading exponents: The fixed point needs to be 
located rather precisely, which can be done with the help of the renormalization group 
trajectories. Bigger lattices should also improve results because then we can get closer to 
the critical surface. 

Our results were obtained with very modest computer resources. All calculations have 
been carried out on a IBM Power PC 250. Once the approximate location of the fixed 
point is found, the MCRG from a 32 3 to a 16 3 system, measuring 14 operators in 10 5 
iterations takes around 50 hours and delivers critical exponents with a ~ 1% precision. 
All the evidence shows that the method is indeed efficient to study the renormalization 
group in the lattice. Our estimates for the critical exponents obtained with quite small 




(5.3) 
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lattices are in perfect agreement with the ones computed by other methods. We think that 
the algorithm is well suited to study more difficult systems (iV-component 4> 4 would be 
interesting) . The application to the crumpling transition in random surfaces is now under 
way. 
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Figure Caption.- 

Fig 1.- Autocorrelation time for the operator (p 2 at A = 6. The Langevin step is At = 0.001. 
The system size is 16 3 . 
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